Metallic phase in the two-dimensional ionic Hubbard model 
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We investigate the phases of the ionic Hubbard model in a two-dimensional square lattice using 
determinant quantum Monte Carlo (DQMC). At half-filling, when the interaction strength or the 
staggered potential dominate we find Mott and band insulators, respectively. When these two 
energies are of the same order we find a metallic region. Charge and magnetic structure factors 
demonstrate the presence of antiferromagnetism only in the Mott region, although the externally 
imposed density modulation is present everywhere in the phase diagram. Away from half-filling, 
other insulating phases are found. Kinetic energy correlations do not give clear signals for the 
existence of a bond-ordered phase. 
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I. INTRODUCTION 

Tight binding models, such as the Hubbard Hamilto- 
nian, have been widely studied for their ability to de- 
scribe different kinds of insulating phases present in con- 
densed matter systems. Interactions between fermions 
give rise to gapped Mott insulators where the fermions 
are immobile at commensurate fillings. For bipartite lat- 
tices, Mott behavior is typically associated with antifcr- 
romagnetic order. On the other hand, the imposition of 
an external periodic potential can drive band insulating 
behavior. In systems where several types of insulating ef- 
fects are present, phase transitions can take place as the 
relative strength of the parameters is changed. Close to 
such transitions, the competition between these different 
terms allows smaller effects, such as kinetic energies, to 
result in low temperature metallic phases or more exotic 
types of correlations such as bond-order. Such compen- 
sating effects are believed to explain enhanced response in 
different quasi-one-dimensional systems such as organic 
materials^, ferroelectric perovskites^ or transition metal- 
oxides^. 

A simple model exhibiting this kind of behavior is the 
ionic Hubbard Hamiltonian 4 , a tight binding model that 
includes an on-site repulsion and a staggered potential 
that takes alternating values on neighboring sites of a 
bipartite lattice. At half-filling, the on-site repulsion 
drives the system into an antifcrromagnetic Mott insu- 
lator (MI), with a homogeneous density. When the stag- 
gered potential dominates, a band insulator (BI) with a 
regular modulation of the density results. A quantum 
phase transition between these two states occurs in the 
phase diagram. 

In one dimension, the ionic Hubbard model has been 
widely studied through a variety of techniques including 
density matrix renormalization group (DMRG), exact di- 
agonalization, bosonization and quantum Monte Carlo 
(QMC) with on-site£&££ or extended^&i&ii^ inter- 
actions. The precise conclusions are still subject to some 
debate, but many of the studies suggest that the direct 
transition from the band insulator to the Mott insula- 



tor is replaced by an intervening insulating bond-ordered 
phase. In two dimensions, dynamical mean field theory 
(DMFT) 13 , cluster-DMFT 14 and determinant quantum 
Monte Carlo (DQMC)^ were used. The existence of 
an intermediate phase between the BI and MI phases 
has been established, as in one dimension, but the pre- 
cise nature of the phase is still under debate (see section 
II). Studies of a bilayer Hubbard modelil and in infinite 
dimension 1 ^ yield similar results with a low temperature 
metal as the intermediate phase. Here the interplanar 
coupling, if sufficiently large, gives a noninteracting band 
structure with a gap at half-filling much like a staggered 
potential. A recent experiment attempts to observe this 
kind of phased. 

A related problem is the study of heterostructures of 
SrTiC>3 and LaTiC>3, one being a MI, the other a BI and 
a metallic phase appearing at the interface between the 
twc3. Some of the phenomena presented here also have 
similarities with those found in corresponding disordered 
models used to describe binary alloys^. 

In this paper, we consider the different phases and 
transitions of the ionic Hubbard Hamiltonian on a two- 
dimensional square lattice using determinant quantum 
Monte Carlo (DQMC), extending our initial studyi^. In 
section II, we introduce the model and summarize our 
previous results. In section III, we briefly describe the 
numerical techniques used while in section IV, we in- 
vestigate further the different phases by measuring the 
charge and magnetic structure factors and the spectral 
function. We then present an analysis of the properties 
away from half filling (section V) and the search for bond- 
order waves (section VI). Conclusions are in section VII. 



II. FERMIONS IN A STAGGERED POTENTIAL 

We consider the ionic Hubbard model 



H = - t ^ (4<r c r' CT + cl, a C ra ) + U ^ n r1 n rl 
(rr')a r 



E 



(n rT +n ri ) 



(1) 



where cJ CT (c r(T ) are the usual fermion creation (destruc- 
tion) operators on site r = (x, y) of a two-dimensional 
square lattice with spin a and where n ra — c^ a c Ta is 
the number operator. The hopping parameter, t, sets 
the energy scale and is fixed at a value of 1, while tun- 
ing the chemical potential \x gives the desired fermion 
filling. The interaction strength U and staggered poten- 
tial (— X)^ x+V ' A/2 are the two competing terms of our 
system, the former encouraging single occupancy of all 
sites and the latter double occupancy of the low energy 
sublatticeJ^ 

Considering the limiting cases U 3> A and A>C/we 
find two distinct insulating phases at half-filling (density 
p = 1). For A = 0, the usual translationally invariant 
Hubbard model arises and, at T = 0, forms an antiferro- 
magnetic Mott insulator at half-filling for any value of U. 
At strong coupling, the gap is set by U, while the effective 
anti- ferromagnetic coupling is set by J = 4£2/U. These 
two energy scales are evident in the magnetic properties 
where T ~ U demarks the onset of moment formation 
on individual sites and T ~ J, the onset of antiferromag- 
netic order between sites. 

In the non interacting limit U = 0, the Hamiltonian 
can be diagonalized in Fourier space (where it reduces to 
a set of independent 2x2 matrices) 



H = J2 Ek c L c ka with E k = ±v/e k 2 + (A/2)2 



k,CT 



where the sum runs over half the Brillouin zone and 
where ek = — 2t (cos(k x ) + cos(k y )) is the energy of a 
free particle of wave vector k = (k x ,k y ). The eigenen- 
ergies E k are then split into two bands separated by a 
gap A and, in this limit, the system is a band insulator 
at half-filling. These two insulating phases can be seen 
schematically in Fig. [TJ 

Recently, it was show n 13 i 14 ' 15 that there is no direct 
transition between these two insulating phases, but that 
an intermediate phase appears (Fig. |2J|. The precise 
nature of the phase is still being debated. In previ- 
ous worki^, we used the determinant quantum Monte 
Carlo algorithm (DQMC) to study the conductivity of 
the model. Performing a finite size and temperature 
analysis, we concluded that the intermediate phase has 
a finite conductivity and is, consequently, metallic. Garg 
et at arrive at the same conclusion in their dynamical 
mean field theory (DMFT) study^ while using a cluster- 
DMFT analysis, Kancharla et al. argue for a bond or- 
dered wave (BOW) phased. 
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FIG. 1: Schematic diagram of electrons in a staggered poten- 
tial. The A sites have lower chemical potential than the B 
sites (the difference being A). For low interaction strength 
U <C A, the electrons prefer to doubly occupy the A sites 
and form a band insulator. For U S> A, it is favorable to 
have only one electron per site and to form a Mott insula- 
tor. The hopping term t then generates an antiferromagnetic 
modulation. 



III. COMPUTATIONAL METHODS 

DQMC? ? is a numerical approach which expresses the 
partition function Z = Tr exp(—pH) as a path integral 
where the inverse temperature (3 is discrctizcd. The in- 
teraction term is decoupled via a Hubbard-Stratonovich 
(HS) transformation, introducing a sum over an Ising 
spin field and leaving the Hamiltonian in a quadratic form 
with respect to the fermion operators. The fermionic de- 
grees of freedom can then be integrated out and the re- 
maining Boltzmann weight is expressed as the determi- 
nant of the product of two matrices that depend on the 
HS spin variables. In the conventional Hubbard model, 
this determinant is positive at half-filling, due to particle- 
hole symmetry. With the staggered potential, the deter- 
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FIG. 2: The phase diagram of the ionic Hubbard mode^. 
Symbols are the result of QMC simulations. The dashed line 
is the strong coupling (t = 0) phase boundary between band 
and Mott insulators. Points (a), (b) and (c) located, respec- 
tively, in the BI, the intermediate phase, and the MI are used 
in Fig. H 
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FIG. 3: Spectral function A(u) for various values of interac- 
tion strength U with A = lona6x6 lattice at = 12 at 
half- filling. The letters (a)-(c) correspond to points on the 
phase diagram Fig. [2] Point (b) exhibits a metallic behavior 
as A(u) is nonvanishing for u = 0. A(ui) exhibits clear gaps 
for points (a) and (c) which are, respectively, band and Mott 
insulating states. 
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FIG. 4: The spectral function A(u>) and conductivity a^c are 
shown as functions of U at /3 — 12 and 16, at half-filling, and 
for A = 0.5. There is a Bl-metal transition at U c i — 0.5 and 
a metal to MI transition at U C 2 — 2.5. 



minant is no longer always positive at half-filling resulting 
in a sign problem which can be severe. 

The remaining sum over the HS field and the evalu- 
ation of observables is achieved by a Metropolis Monte 
Carlo algorithm. It is difficult to obtain good estimates 
for the observables when the determinant is often nega- 
tive, i.e., when the mean value of the sign of the determi- 
nant is small. This problem is important for intermediate 
values of U and A, while it is not severe in the limiting 
cases where we retrieve simple BI or MI phases (U or A 
small or large). 

With this formulation, space and (imaginary-)time 
separated Green functions 

G(R,r) = ( Cr+Rff (r) C t (T (0)) 



are easily expressed as elements of the matrices. ^From 
the spatial Fourier transform of G(R, r), one can extract 
the spectral function (essentially the density of states), 
by inverting the Laplace transform 



G(k = 0,r) 



1 



(2) 



Solving for A(uj) remains a non-trivial task that is per- 
formed with a method of analytic continuation 23 . We 
measure correlation functions for charge and spin den- 
sity waves (CDW and SDW, respectively), given by 



Ccd\ 
sd\ 



.(R) 



-Y, 



([n r+ Rt ± »v+m] K-T ± n rl \) 



where + corresponds to CDW and — to SDW. The stag- 
gered potential imposes a checkerboard modulation on 
the density, not only in the BI phase but in all phases 
as it breaks explicitly the translational invariance. An- 
tiferromagnetic magnetic correlations are also expected 
in the MI phase. Thus, for both types of correlations, 
we expect ordering with a wave vector k w = (w, tt), char- 
acterized by the divergence of corresponding structure 
factors 



S cdw 
sdw 



R 



R 



G c dw (R) 

sdw 



(3) 



IV. 



SPECTRAL FUNCTION, MAGNETIC AND 
CHARGE ORDERS 



To complete our previous analysis^ of the intermedi- 
ate phase shown in Fig. [2] we evaluated the charge and 
spin correlations and spectral function. When the size is 
not explicitly mentioned, we used aiV = 6x6 = 36 sites 
lattice. N x = yN denotes the side of the lattice. 

Fig. [3] shows the spectral function for A = 1 (along 
the vertical line shown in Fig. \2§ at (3 = 12 and for three 
values of U corresponding to the three different phases. 
For U = 0.5 and [7 = 4, corresponding to the BI and MI 
phases, the spectral function shows a gap around u> = 0, 
thus indicating an insulator. In contrast, the spectral 
function for U — 2 is non zero and shows a quasi-particle 
peak at lu — which indicates metallic behavior. 

These results are consistent with measurements of the 
conductivity in our previous article, as is confirmed by 
the direct comparison of the spectral function at zero 
frequency A{ui = 0) and the conductivity a^ c shown in 
Fig. [4] For A = 0.5 and as U is increased, the behaviors 
of these two quantities are similar: they both show, for 
intermediate values of U , a peak that grows with (3. The 
ranges of interaction where there are peaks are slightly 
different because of finite size effects. 

U c i and U C 2 are the critical repulsion energies that 
cause BI to metal and metal to MI phase transitions for a 
fixed value of A. These critical values are obtained by a 
finite temperature analysis of the conductivity 15 . As the 
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FIG. 5: Top: Antiferromagnetic structure factor S s dw as a 
function of repulsion U shows antiferromagnetic order devel- 
oping for large U. Bottom: Conversely, the charge structure 
factor Scdw decreases as U increases. A = 0.25 in both figures. 



temperature is lowered, the conductivity increases in a 
metal and decreases in an insulator. Locating the points 
where the behavior of the conductivity changes gives the 
critical values. For A = 0.5, the Bl-metal transition is 
located at U c \/t ~ 0.5, the metal-MI transition is located 
at U C 2/t ~ 2.5 and the metallic region is in between. 

Structure factors S s d w and S c d w are shown in Fig. [5] 
for A= 0.25 and for inverse temperatures j3 = 10, 12 and 
16. For these values of A, the metallic boundaries are 
U c i ~ 0.1 and U C 2 — 1.7. As U is increased, S s d w in- 
creases as the system is driven towards the MI phase, 
while S'cdw decreases, in moving away from the BI phase. 
The transitions in these quantities when we cross the 
phase boundaries are not sharp on lattices of fixed, finite 
size owing to finite size rounding effects. In addition, as 
stated above, a smoother evolution is expected for S'cdw 
as a modulation of p is imposed by the external potential 
for all parameter values. 

In a conventional Hubbard model, the MI at half fill- 
ing is always an antiferromagnet. However, it has been 
suggested recently for related models^ that the metal- 
MI and para- to antiferro-magnetic (PM-AF) transitions 
may not take place together. A finite size scaling anal- 
ysis of S'sdw for different values of U is needed to check 



0.2 



0.15 



0.1 



CO 



0.(15 



o 


u= 


=0.25 


□ 


u= 


=1.50 


o 


u= 


=2.25 


A 


u- 


=2.50 


< 


u- 


=2.75 


V 


v= 


=3.00 



A- 



0.05 



8 8 



1 



-H" 
-■©- 



□ 

o 



0.1 

1/7Y 



0.15 



FIG. 6: The antiferromagnetic structure factor S a dvi/N as a 
function of inverse lattice size 1/N X at A = 0.5 and /3 = 
12 for various U. The different values of U correspond to 
different phases: U = 0.25 to the BI, U = 1.50, 2.25 and 2.5 
to the metal, and U = 2.75 and 3 to the MI. The lines are 
least-squares fit to the data and the extrapolation leads to 
the value of S s d w = m2/3 for an infinite system. There is 
antiferromagnetic long range order only in the MI phase. 



whether these two boundaries coincide in our model. The 
structure factor S s d w increases like the number of sites N 
in the thermodynamic limit 



^sdw 



m 

T 



o 



i 



(4) 



where m2 is the antiferromagnetic order parameter. 
Fig. [6] shows structure factors for sizes up to N = 12x12 
sites for different values of U. A is set to 0.5 and the 
metallic phase appears between U c i — 0.5 and U C 2 — 2.5. 
We apply a least squares fit to the data, and extrapolate 
to an infinite system. We find a non-vanishing structure 
factor only in the MI phase (U = 2.75 and 3). In all the 
other cases, including points located close to the metal- 
MI transition (U — 2.25 and 2.5), the structure factor 
goes to zero as the size is increased. 

The PM-AF transition, therefore, takes place between 
U = 2.5 and U — 2.75 and appears to be simultaneous 
with the metal-MI transition. Although our study on 
limited sizes cannot completely exclude completely the 
possibility of an exotic AF metallic phase located around 
the metal-MI boundary, our data suggest that it does not 
occur. 



V. AWAY FROM HALF-FILLING 

We extend our study to determine the behavior away 
from half-filled band. In Fig. [71 the density p is plotted 
as a function of chemical potential. The compressibility 
k is given by dp/ dp. For k—0, incompressible insulating 
phases are present which is indicated by plateaux in the 
p versus p curve. Both the MI (U = 2.0, A = 0.2) and BI 
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FIG. 7: Density, p, as a function of chemical potential p for 
the different phases. The Mott (U = 2.0, A = 0.2, circles) 
and band (U — 0.5, A = 1.0, diamonds) insulating phases 
have plateaux as the compressibility k = dp/dp = at half- 
filling. For the metallic phase (U = 1.0, A = 0.5, squares), 
there is no plateau in the region around p = 1. 




FIG. 8: Spin (top) and charge (bottom) structure factors as 
functions of density. The parameters are the same as in Fig. [7] 
As p approaches half-filling, S'sdw reaches its maximum value 
for the Mott phase. Conversely, S c dw is maximum at p = 1 
for the band insulator. The metallic phase has a slight peak 
in S'cdw but not in S s d„. 



(U = 0.5, A = 1.0) phases have this distinctive plateau 



at half-filling. The metallic nature of the intermediate 
phase is once again confirmed by the absence of a plateau 
around p = 1. 

Using the same parameters U and A, we show in Fig. [5] 
the evolution of the structure factors as functions of p. 
As expected, for the parameters corresponding to the MI 
at p = 1, S'sdw shows a strong peak around half-filling 
indicating the presence of the AF-MI in this case and 
its rapid suppression upon doping. In this case, charge 
correlations always remain small, as double occupancy of 
sites is essentially forbidden. Similarly, for the parame- 
ters corresponding to the BI, the peak is present in S c d w 
at p = 1 and S'sdw remains small for all densities, as the 
occupation of larger energy site is suppressed. For the 
metallic case, there is no peak in the magnetic structure 
factor at p = 1 but a moderate peak in S c d w - 

In all these three cases, away from half-filling, the 
structure factors are rapidly suppressed and the system 
is compressible and becomes a paramagnetic metal. 
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FIG. 9: The density, p, as a function of chemical potential 
p for A = 8 and U = 4 at f3 = 10. The plateaux represent 
the various insulating phases. The regions where ft > are 
metallic phases. 

For stronger interactions and staggered potential, 
other ordered phases can be found at quarter and three 
quarters fillings. For U = 4 and A = 8, there is a band 
insulator at half-filling shown by the p — 1 plateau in 
Fig. [5] Additional plateaux emerge at p = 0.5 and 1.5. 
The plateau at p = 0.5 corresponds to a state where 
all the lower energy sites are occupied by one fermion. 
Adding a fermion to an already occupied site increases 
the energy approximately by U, while adding it to an 
empty (high energy) site increases it by A. There is then 
an energy gap, roughly equal to the smaller of the two 
parameters (U in this case). We then have a Mott in- 
sulator at quarter filling and, symetrically, another one 
at three quarters filling. In the case where U is larger 
than A, the situation is reversed and we expect a MI at 
half-filling and a BI for p — 0.5 and p — 1.5 (this case 
is difficult to treat numerically because of a strong sign 
problem) . 

Although the system is incompressible for p = 0.5, 
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FIG. 10: (a) Antiferromagnetic structure factor S s d v , (b) 
charge structure factor S c< j w and (c) conductivity, <7d c as func- 
tions of p at U — 2, A = 8 and /3 = 10. Note that Ode vanishes 
in the three solid phases, p = 1/2, 1, 3/2. 



there is no concommitant magnetic order: in Fig. I10[ 
S'sdw always remains small as p is varied away from half- 
filling and there is no AF peak around quarter filling 
(the maximum being around p = 0.75). The AF or- 
dering appears when the fermions jump between neigh- 
boring sites, thus creating an effective antiferromagnetic 
coupling between spins. At quarter filling, the particles 
are surrounded by empty sites and there is no effective 
coupling. As for S c dw we only observe a peak at half- 
filling (Fig. HUJl. A plot of the conductivity (Fig. Hol e)) 



versus p shows that, indeed, <7dc vanishes for p = 1/4 
and p = 3/4 corresponding to the plateaux observed in 
Fig. m 

VI. BOND-ORDERED WAVES 

Long range Bond-ordered waves (BOW) are known to 
exist in one dimension in the ionic Hubbard model^^ as 
well as in other one-dimensional systems exhibiting den- 
sity modulation such as extended Hubbard models 10,11 . 
(In the latter case, the two numerical treatments, us- 
ing the stochastic series expansion** and density matrix 
renormalization group* -, disagree on the precise width of 
the BOW phase as well as over what range of coupling- 
values it is present.) It was suggested by Kancharla et 
al)^ that the intermediate phase observed at half-filling 
in the two-dimensional case studied here was in fact a 
bond-ordered state, as in one dimension. Other stud- 
ies have attempted to observe bond-order in two dimen- 
sions, for example in extended Hubbard models*^. In 
two dimensions, bond orders have mostly been obtained 
in Hamiltonians possessing rather special geometries or 
interactions, e.g. on checkerboard lattices*^ or by includ- 
ing additional terms such as ring exchanges*^. 
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FIG. 11: Schematic diagram of how a bond-ordered wave oc- 
curs in one dimension. For U ~ A both large, the electrons 
may move to a site of higher chemical potential as shown 
with the spin-up electron. The spin-down electron in the left 
well is not likely to move to the right because it will cost a 
large energy A to doubly occupy the site of higher chemical 
potential. Prohibited movement to left causes a decrease in 
kinetic energy for that bond. This alternation of two bonds 
with weak and strong kinetic energy bonds is an indication of 
the appearance of a bond order. Strong bonds are continuous 
thick lines, weak bonds are dotted lines. 

Hence, bond order seems in general more difficult to 
obtain in two dimensions than in one. In our system, 
we can understand this with the argument sketched in 
Fig. QT] Starting from a situation where all the particles 
are located on the low energy sites (A sites), a jump to 
a high energy site (B site) increases the energy roughly 
by A — U . If A and U are of the same order of magni- 
tude, this energy cost is small and the jump has a high 
probability to happen. This B site is connected to two 
A sites: the one where the jumping particle comes from 
(to the right in Fig.fTTj) and another one (to the left) still 
occupied by two particles. The jumps on the left bond 
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FIG. 12: The bond-ordered susceptibility as a function of U 
at p = 1 for AT = 12 x 12 and /3 = 10. The vertical lines 
separate the three different phases for two different values of 
the staggered potential. The maximum of this bond-order 
parameter is in the metallic region but is of order 10 -4 . 
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FIG. 13: The bond-ordered susceptibility and conductivity as 
functions of p for A = 1, U — 0.5, and /3 = 12. Doping away 
from half-filling, the susceptibility is non-zero even when the 
band is only partially full and the system is unquestionable 
metallic. 



are suppressed: the particle on the B site and the par- 
ticle having the same spin (an up spin in the sketch) on 
the A site cannot jump because of Pauli principle, while 
a jump of the remaining particle has a large energy cost 
A. 

In short, having a particle jumping on a bond (and 
thus a contribution to the kinetic energy on this bond) 
will inhibit jumps on the following bond. In one dimen- 
sion, this can lead to a bond order with alternate strong 
and weak values of the kinetic energy. On a two dimen- 
sional square lattice, the appearance of a strong bond 
will inhibit the three other bonds linked to the same site 
and it seems more difficult to generate regular arrays of 
weak and strong bonds. 



To check if there is a BOW, we examined the correla- 
tion function of the kinetic energy 

r 

where k% is the kinetic energy operator connecting sites 
r and r + x (x being a unit vector directed along the x 
direction), 

k r = y ] ( C r+x,er Cr ,°' c r.<y Cr +x,<y)- 

The Fourier transform of K XX (R) always exhibits a very 
small peak for = (it, tt), in all the different phases 
at or away from half-filling. Finite size scaling analysis 
indicates that these peaks persist in large size lattices, 
but this is the case for all the different phases, not only 
in the intermediate phase. We also studied the BOW 
susceptibility 11 for this wave vector 

Xbow= £ e lk - R (fc? +R (r)^(0)) (5) 
R,r,T 

Fig. [12] shows the evolution of Xbow at half-filling as a 
function of U and for two different values of A. The 
curve is maximum in the intermediate phase but is also 
non zero in other phases. Similar results are obtained for 
bonds oriented along the y direction and for mixed x-y 
cases. 

We do not believe these signals to indicate the pres- 
ence of a BOW phase for the following reasons (a) The 
values obtained are extremely small, (b) These "BOW 
signals" are ubiquitous in all three phases and not spe- 
cific to the intermediate one. The susceptibility is, for 
example, nearly as important in the MI phase as in the 
intermediate phase. As soon as U and A are both non- 
zero, we observe this very weak non-zero value of Xbow 
(c) Starting from the metallic phase (U — 0.5, A = 1) 
and doping away from half-filling the system is certainly 
known to be metallic (see Fig. [T5|) yet the BOW suscep- 
tibility is non zero and comparable to the value found at 
half-filling. For all these reasons, we believe that these 
tiny signals are only a weak signature of an interplay 
between the interaction and externally imposed checker- 
board potential. 

VII. CONCLUSIONS 

We have extended our previous study of the phase di- 
agram of the ionic Hubbard model 15 . In addition to new 
QMC results for the conductivity, we calculated the spec- 
tral function, the kinetic energy correlation function (the 
order parameter for the BOW phase) and the magnetic 
and charge order parameters. Our results point to the 
presence of a normal conducting metallic phase interven- 
ing between the band and Mott insulating phases. In 
particular, we have not found a clear signal in the BOW 
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order parameter and susceptibilities indicating the pres- 
ence of such a phase. 

We also studied the system doped away from half fill- 
ing. When the band or Mott insulating phase is doped, a 
conducting phase is obtained except at fillings of 1 /4 and 
3/4 where we found insulating phases without magnetic 



order. 
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